/* Examples/jerasure_03.cpp
 * James S. Plank

Jerasure - A C/C++ Library for a Variety of Reed-Solomon and RAID-6 Erasure Coding Techniques
Copright (C) 2007 James S. Plank

James S. Plank
Department of Electrical Engineering and Computer Science
University of Tennessee 
Knoxville, TN 37996
plank@cs.utk.edu
*/

/*
 * $Revision: 1.2 $
 * $Date: 2008/08/19 17:41:40 $
 */
    

#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "jerasure-2.h"

#define talloc(type, num) (type *) malloc(sizeof(type)*(num))

void usage(char *s)
{
  fprintf(stderr, "usage: jerasure_03 k w - Creates a kxk Cauchy matrix in GF(2^w). \n\n");
  fprintf(stderr, "       k must be < 2^w.  Element i,j is 1/(i+(2^w-j-1)).  (If that is\n");
  fprintf(stderr, "       1/0, then it sets it to zero).  \n");
  fprintf(stderr, "       It then tests whether that matrix is invertible.\n");
  fprintf(stderr, "       If it is invertible, then it prints out the inverse.\n");
  fprintf(stderr, "       Finally, it prints the product of the matrix and its inverse.\n");
  fprintf(stderr, "       \n");
  fprintf(stderr, "This demonstrates: JER_Matrix::Print()\n");
  fprintf(stderr, "                   JER_Matrix::Invertible()\n");
  fprintf(stderr, "                   JER_Matrix_Inverse()\n");
  //fprintf(stderr, "                   jerasure_matrix_multiply().\n");
  if (s != NULL) fprintf(stderr, "%s\n", s);
  exit(1);
}

int main(int argc, char **argv)
{
  unsigned int k, w, i, j, n;
	JER_Matrix *jm;
	JER_Matrix *identity;

  if (argc != 3) usage(NULL);
  if (sscanf(argv[1], "%d", &k) == 0 || k <= 0) usage((char *)"Bad k");
  if (sscanf(argv[2], "%d", &w) == 0 || w <= 0 || w > 31) usage((char *)"Bad w");
  if (k >= (1 << w)) usage((char *)"K too big");

	jm = new JER_Matrix(k, k, w);

  for (i = 0; i < k; i++) {
    for (j = 0; j < k; j++) {
      n = i ^ ((1 << w)-1-j);
			jm->Elts[i][j] = (n == 0) ? 0: galois_single_divide(1, n, w);
    }
  }

  printf("The Cauchy Matrix:\n");
	jm->Print();
  i = jm->Invertible();
  printf("\nInvertible: %s\n", (i == 1) ? "Yes" : "No");
  if (i == 1) {
    printf("\nInverse:\n");
		i = jm->Inverse();
		jm->Print();

		identity = new JER_Matrix(k, k, w);

		identity->CreateIdentity();
    printf("\nInverse times matrix (should be identity):\n");
		identity->Print();
  }
  return 0;
}

